function N = reference_basis_function_2d(basis_type)
switch basis_type
    case "P1"
        N = cell(3,1);
        N{1} = @(x,y) 1 - x - y;
        N{2} = @(x,y) x;
        N{3} = @(x,y) y;
    case "P1b"
        N = cell(4,1);
        N{1} = @(x,y) 1 - x - y;
        N{2} = @(x,y) x;
        N{3} = @(x,y) y;
        N{4} = @(x,y) -x.*y.*(27*x + 27*y - 27);
    case "P2"
        N = cell(6,1);
        N{1} = @(x,y) 2*x.^2 + 4*x.*y - 3*x + 2*y.^2 - 3*y + 1;
        N{4} = @(x,y) -4*x.*(x + y - 1);
        N{2} = @(x,y) x.*(2*x - 1);
        N{5} = @(x,y) 4*x.*y;
        N{3} = @(x,y) y.*(2*y - 1);
        N{6} = @(x,y) -4*y.*(x + y - 1);
    otherwise
        error("Invalid basis type.");
end
end